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Methods for identification of dynamical patterns in networks suffer from effects of arbitrary time 
scales that need to be imposed a priori. Here we develop a principled method to identify patterns 
on dynamics that take place on network systems, as well as on the dynamics that shape the network 
themselves, without requiring the stipulation of relevant time scales, which instead are determined 
solely from data. Our approach is based on a variable-order hidden Markov chain model that 
generalizes the stochastic block model for discrete time-series as well as temporal networks, without 
requiring the aggregation of events into discrete intervals. We formulate an efficient nonparametric 
Bayesian framework that can infer the most appropriate Markov order and number of communities, 
based solely on statistical evidence and without overfitting. 


I. INTRODUCTION 

To reveal the mechanisms that form complex systems, 
community-detection methods describe large-scale pat¬ 
terns in their networks of interactions [T]. Tradition¬ 
ally, these methods describe only static network struc¬ 
tures, without taking into account that networks form 
dynamically over time. Only recently have researchers 
proposed methods that incorporate higher-order tempo¬ 
ral structures to capture dynamics of networks BHiD]. 
Another class of methods attempt to describe the dy¬ 
namics that take place on a network, which can also be 
used as a proxy for the network structure. Traditionally, 
these methods use memory less Markov chains |21H23| . 
but researchers have recently showed that incorporating 
memory effects can reveal crucial information about the 
system |24fE5] . However, both avenues of research have 
encountered central limitations: On the one hand, meth¬ 
ods that uses memory to describe dynamics on networks 
rely on extrinsic methods to detect the appropriate mem¬ 
ory order On the other hand, methods that 

attempt to describe dynamics of networks adapt static 
descriptions by aggregating time windows into discrete 
layers |28| . and completely ignore dynamics within the 
time windows. 

Both for dynamics on and of networks, the methods 
require or impose ad hoc time scales that should instead 
be determined solely from data, or avoided entirely. Fur¬ 
thermore, since most descriptions also depend on a large 
number of additional degrees of freedom from all pos¬ 
sible network partitions, they are prone to overfitting 
when random fluctuations in high dimensional data are 
mistaken for actual structure |29| . The imposed ad hoc 
time scales exacerbate this problem, since they further 
increase the number of degrees of freedom in the descrip- 
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tions. Without a principled methodology to counteract 
this increasing complexity, it becomes difficult to sepa¬ 
rate meaningful effects from artifacts. 

We present a general and principled approach to this 
problem by tackling dynamics on and of networks simul¬ 
taneously (see Fig. [^. In contrast to approaches that in¬ 
corporate temporal layers in methods for static network 
descriptions, we build our approach on describing the ac¬ 
tual dynamics. We first formulate a generative model 
of discrete temporal processes based on arbitrary-order 
hidden Markov chains with community structure |30f[55] . 
Since our model generates arbitrary sequences of events, 
it does not require aggregation in time windows, and 
hence needs no a priori imposition of time scales. This 
model can be used to describe dynamics taking place 
on network systems that take into account memory ef¬ 
fects BUBS] of arbitrary order. We then use this model 
as the basis for the description of temporal networks, 
where the sequence of events represents the occurrence 
of edges in the network [8]. In both cases, we employ a 
nonparametric Bayesian inference framework that allows 
us to select the most parsimonious model according to the 
statistical evidence available, and hence is able to detect 
the most appropriate Markov order, as well as other prop¬ 
erties such as the number of communities, while avoiding 
overfitting. In particular, if there is no structure in the 
data — either via a fully random dynamics or a lack of 
large-scale structure in the network — our method is able 
to identify this. As we also show, the model can be used 
to predict future network dynamics and evolution from 
past observations. 

Our model builds on the original idea of the stochastic 
block model (SBM) BHISSIj where the nodes are divided 
into groups, and the edge placement probabilities depend 
on the node memberships. In Sec. |TTjwe adapt this idea 
to discrete-time Markov chains, that operate on condi¬ 
tional probabilities in sequences of tokens, where both 
the individual tokens and their preceding subsequence 
are divided into groups. Then, in Sec. |III[ we extend this 



2 


(a) Dynamics on networks 



Temporal sequence of tokens (nodes in a network): 
{xt} = "Ituwasutheubestuofutimes" 


(b) Dynamics of networks 

Edge 


1 I I i I 

2 3 4 5 6 7 Tiine 

Temporal sequence of edges in a network: 

{xt} = {(1, 2), (4, 3), (5, 2), (10, 8), (7, 2), (9, 3), (3, 4)} 

Figure 1. Our modeling framework allows simultaneously for 
the description of (a) arbitrary dynamics taking place on net¬ 
works, represented as sequence of arbitrary “tokens” that are 
associated with nodes, and (b) dynamics of networks them¬ 
selves, where the tokens are edges (pairs of nodes) that appear 
in sequence. 
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model to a community-based temporal networks that op¬ 
erates on the sequence of added edges in the evolution of a 
network. In both cases, we illustrate the use of the model 
with empirical data. In Sec. |IV| we discuss extensions of 
the approach to dynamical systems with continuous time, 
as well as nonstationary Markov chains. We finalize in 
Sec. El with a conclusion. 


II. INFERENCE OF MARKOV CHAINS 


Here we consider general time-series composed of a se¬ 
quence of discrete tokens {a;t}, where Xt is a single token 
from an alphabet of size N observed at discrete time t, 
and Xt-i = {xt-i, ■ ■ ■ ,Xt-n) are the previous n tokens 
at time t. An nth-order Markov chain with transition 
probabilities p{xt\xt-i) generates such a sequence with 


probability 

Pi{xt}\p) = J|p(a:t|5't_i), = J|p(a;|f)“-'®, (1) 

t X,X 

where is the number of transitions x ^ x in {a^t}. 
Given a specific sequence {xt}, we want to infer the tran¬ 
sitions probabilities p{x\x). The simplest approach is to 
compute the maximum-likelihood estimate, i.e. 

pix\x) = argmaxP({a:t}|p) = (2) 

p{x\x) ^x 


where as = ax,x, which amounts simply to the fre¬ 
quency of observed transitions. Putting this back into 
the likelihood of Eq. we have 

lnP({a:t}|p) = ^aa;,5rln—■ (3) 

_ ^X 

x,x 

This can be expressed via the conditional en- 
tropy H(XIX) = - J^sf(x) jd(xlx) Inp(xlx), with 
lnP{{xt}\p) = -EH{X\X), where E = J2x,x^x,s is the 
total number of observed transitions. Hence, the max¬ 
imization of the likelihood of Eq. yields, seemingly, 
the transition probabilities which most compress the se¬ 
quence. There is, however, an important caveat with this 
approach. Namely, it cannot be used when we are inter¬ 
ested in determining the most appropriate order n of the 
model. This is because if we increase n, the maximum 
likelihood of Eq. will also invariably increase: As the 
number of memories becomes larger, while the total num¬ 
ber of transition remains fixed, the conditional entropy 
estimated from the frequency of transitions is bound to 
decrease. Hence, for some large enough value of n there 
will be only one observed transition conditioned on every 
memory, yielding a zero conditional entropy and a max¬ 
imum likelihood of one. This would be an extreme case 
of overfitting, where by increasing the number of degrees 
of freedom of the model one is not able to distinguish 
actual structure from stochastic fluctuations. Also, this 
approach does not yield true compression of the data, 
since it omits the description of the model (which be¬ 
comes larger with n), and thus is crucially incomplete. 
To address this problem, one must use a Bayesian for¬ 
mulation, and maximize instead the complete evidence 

Pi{xt}) = J dp P{{xt}\p)P{p), (4) 

which consists in the sum of all possible models weighted 
according to some prior probabilities P{p), that encode 
our a priori assumptions. This approach has been con¬ 
sidered before by Strelioff et al [33], where it has been 
shown to yield the correct model order for data sampled 
from Markov chains — as long as there is enough statis¬ 
tics — as well as meaningful values also when this is not 
the case, that balances the structure present in the data 
with its statistical weight. 
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Figure 2. Schematic representation of the Markov model 
with communities for a sequence of tokens {xt} = 
"Ituwasutheubestuofutimes". The nodes are either memo¬ 
ries (top row) or tokens (bottom row) and the directed edges 
represent observed transitions, (a) A partition of the tokens 
and memories for a n = 1 model, (b) A “unified” formulation 
of a n = 1 model, where the tokens and memories have the 
same partition, and hence can be represented as a single set 
of nodes, (c) A partition of the tokens and memories for a 
n = 2 model. 

Although this Bayesian approach addresses satisfacto¬ 
rily the overfitting problem, it misses opportunities of 
detecting structnres in data. As we show below, it is 
possible to extend this model in snch a way as to make a 
direct connection to the problem of finding commnnities 
in networks, yielding a stronger explanatory power when 
modeling sequences, and serving as a basis for a model 
where the sequence itself represents a temporal network. 

A. Markov chains with communities 

Instead of directly inferring the transition probabilities 
of Eq. ^ we propose an alternative formnlation: We as- 
snme that both the memories and tokens are distributed 
in disjoint groups (see Fig. [^. That is, bx G 
and [B]\[ + 1, Bn + Bm] are the group memberships 
of the tokens and memories, respectively, snch that the 
transition probabilities can be parametrized as 

p{x\x) = OxXb^b^- (5) 

Here 9x is the relative probability at which token x is se¬ 
lected among those that belong to same gronp (playing a 
role analogous to degree-correction in the SBM m) m, 
and Xrs is the overall transition probability from mem¬ 
ory gronp s to token group r. In the case n = 1, for 
example, each token appears twice in the model, both 
as token and memory. (An alternative and often useful 
approach for n = 1 is to consider a single unified parti¬ 
tion for both tokens and memories, as shown in Fig. 
and described in detail in Appendix]^) The maximum 
likelihood estimates for the parameters are 

= ex = —, ( 6 ) 

where e^s is the number of observed transitions between 
groups r and s, ^rs is the total outgoing (or 


Previous n = 3 airports, x 



Atlan^ Las Vegas airport, O; 


Figure 3. Part of the US Air flights itineraries dur¬ 
ing 2011. Only itineraries containing Atlanta or Las Ve¬ 
gas are shown. Edges incident on memories of the type 
X = (it-i, Atlanta,xt-a) are colored in red, whereas x = 
(xt-i,Las Vegas,xt-a) are colored in blue. The node colors 
and overlaid hierarchical division correspond to part of the 
n = 3 model inferred for the whole dataset. 

incoming) transitions, and kx is the total number of oc¬ 
currences of token x. Pntting this back in the likelihood, 
we have 

lnP({a;t}|6, In + '^kxlnkx. (7) 

r<s X 

This is almost the same as the maximum likelihood of the 
degree-corrected stochastic block model (DCSBM) |35| . 
where ax^x plays the role of the adjacency matrix of a bi¬ 
partite multigraph connecting tokens and memories. The 
only differences are constant terms that do not alter the 
position of the maximum with respect to the node parti¬ 
tion. This implies that, in certain sitnations, there is no 
difference between inferring the structure directly from 
its topology or from dynamical processes taking place on 
it (see Appendix [P] for more details). 

As before, this maximnm likelihood approach cannot 
be used if we do not known the order of the Markov 
chain, otherwise it will overfit. In fact, this problem 
is now aggravated by the larger number of parameters 
this model has. Therefore, we employ a Bayesian for¬ 
mnlation and construct a generative processes for the 
model parameters themselves. We do this by introducing 
prior probability densities for the parameters 'Dr{{9x}\a) 
and I?s({Ar.s}|/3), with hyperparameter sets a and /3, and 
computing the integrated likelihood 

P{{xt}\a,l3,b)= J d0dA Pi{xt}\b,X,e) 

xl[Vr{{Ox}\a)l[Vs{{Xrsm). ( 8 ) 
r s 

Now, instead of inferring the hyperparameters, we can 
make a noninformative choice for a and /? that reflects 
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our a priori lack of preference towards any particular 
model [3H] • Doing so in this case yields a likelihood (see 
Appendix [A] for details), 


P{{xt}\b, {ej) = P{{xt}\b, {crs}, {k;^}) 


where 




P{{xt}\b,{ers},{k.}) = 

P{{kx}\{ers},b) = 

^({e..}|{e.}) = n 


Ur eriJls 

n 


B 


N 


6rs}|{Cs}), 

(9) 

tHw, 

(10) 

1 

(11) 


(12) 


with ((™)) = expression above 

has the following combinatorical interpretation: 
P{{xt}\b,{ers},{kx}) corresponds to the likelihood 
of a microcanonical model |39| where a random sequence 
{xt} is produced with exactly Crs total transitions 
between groups r and s, and with each token x occurring 
exactly times (see the Appendix for a proof). 
The remaining likelihoods are the prior probabilities 
on the discrete parameters {crs} and {kx}, which are 
uniform distributions of the type 1/fl, where D is 
the total number of possibilities given the imposed 
constraints m- 

If we apply Stirling’s factorial approximation to Eq. 
we obtain \n P{{xt}\b, {crs}, {kx}) « \nP{{xt}\b,X,6), 
with the right-hand side given by Eq. Hence, the re¬ 
maining terms (Eqs. 11 and[T^ serve as a penalty to the 
maximum likelihood estimate that prevents overfitting as 
the size of the model increases via n, Bn, or Bm |43H46| . 

To make the above model fully nonparametric, we in¬ 
clude priors for the remaining discrete parameters: the 
node partitions {bx} and {b^}, as well as the memory 
group counts, {cg}, so that the complete joint likelihood 
becomes 


P{{xt}, b, {ej) = P{{xt}\b, {eg})P{{b,})P{{bs})P{{es}) 

(13) 

This can be done in the same manner as for the SBM |43h 
145] . In particular, to avoid underfitting the model |43| 
and to reveal hierarchical modular structures [H], the 
uniform priors of Eqs. El and can be replaced by 
multilevel Bayesian hierarchies (i.e. a sequence of pri¬ 
ors and hyperpriors). In order to fit the model above 
we need to find the partitions {bx} and {bg} that max¬ 
imize the joint likelihood P({xt},b), or, equivalently, 
minimize the description length E = —hiP{{xt},b) = 
— lnP({xt}|5) —lnP(&) I351I37]- This latter quantity has 
a straightforward but important information-theoretical 
interpretation: it corresponds to the amount of informa¬ 
tion necessary to describe both the data and the model 
simultaneously. Because of its complete character, mini¬ 
mizing it indeed amounts to achieving true compression 


of data, differently from the parametric maximum like¬ 
lihood approach mentioned earlier. Because the whole 
process is functionally equivalent to inferring the SBM 
for networks, the same algorithms can be used |3H] (see 
Appendix [A| for a summary of the inference details). 

This Markov chain model with communities succeeds 
in providing a better description for a variety of empiri¬ 
cal sequences when compared with the common Markov 
chain parametrization (see Table |^. Not only do we ob¬ 
serve a smaller description length systematically, but we 
also find evidence for higher order memory in all exam¬ 
ples. We emphasize that we are protected against over¬ 
fitting: If we randomly shuffle the order of the tokens 
in each dataset, we always infer a fully random model 
with n = 1 and Bn = Bm = 1. Therefore, we are not 
susceptible to the spurious results of nonstatistical meth¬ 
ods [^ . 

To illustrate the effects of community structure on the 
Markov dynamics, we use the US air flight itineraries as 
an example (Fig. |^. In this dataset, the itineraries of 
1,272,696 passengers were recorded, and here each air¬ 
port stop is treated as a token in a sequence (see Ap¬ 
pendix for more details). When we infer our model, 
the itinerary memories are grouped together if their des¬ 
tination probabilities are similar. As a result, it becomes 
possible, for example, to distinguish “transit hubs” from 
“destination hubs” [23|. We use Atlanta and Las Vegas 
to illustrate: Many roundtrip routes transit through At¬ 
lanta from the origin to the final destination and return 
to it two legs later on the way back to the origin. On the 
other hand. Las Vegas often is the final destination of 
a roundtrip such that the stop two legs later represents 
a more diverse set of origins (Fig. |^. This pattern is 
captured in our model by the larger number of memory 
groups that involve Las Vegas than those that involve 
Atlanta. It is indeed these equivalence classes (i.e. mem¬ 
ories and tokens that belong to the same group) that 
effectively reduce the overall complexity of the data, and 
provide better compression at a higher memory order, 
when compared to the traditional Markov chain formula¬ 
tion. Consequently, the community-based Markov model 
can capture patterns of higher-order memory that con¬ 
ventional methods obscure. 

As we show in the next section, this model can also be 
used as the basis for a temporal network model, where 
the tokens in the sequence are edges placed in time. 
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Table I. Description length E = — log 2 P({a:t}, &) (in bits) as well as inferred number of token and memory groups, Bm and 
Bm, respectively, for different data sets and different Markov order n (see Appendixfor dataset descriptions). The value 
E' = — logj P({ 2 ;(}) corresponds to the direct Bayesian parametrization of Markov chains of Ref. |36| . Values in grey correspond 
to the minimum of each column. At the bottom is shown the compression obtained with gzip and LZMA, two popular variations 
of Lempel-Ziv [411142| . for the same datasets. 


III. TEMPORAL NETWORKS 


A general model for temporal networks consists in 
treating the edge sequence as a time series ElEllin]. We 
can in principle use the model above without any mod¬ 
ification by considering the observed edges as tokens in 
the Markov chain, i.e., Xt = where i and j are the 

endpoints of the edge at time t (see Fig. &)• However, 
this can be suboptimal if the networks are sparse, i.e., if 
only some relatively small subset of all possible edges oc¬ 
cur, and thus there are insufficient data to reliably fit the 
model. Therefore, we adapt the model above by including 
an additional generative layer between the Markov chain 
and the observed edges. We do so by partitioning the 
nodes of the network into groups, i.e. Ci € [IjC] deter¬ 
mines the membership of node i in one of C groups, such 
that each edge (i,J) is associated with a label (ci,cj). 
The idea then is to define a Markov chain for the se¬ 
quence of edge labels, with the actual edges being sam¬ 
pled conditioned only on the labels. Since this reduces 
the number of possible tokens from 0{N‘^) to 0{C^), it 
has a more controllable number of parameters that can 
better match the sparsity of the data. We further assume 
that, given the node partitions, the edges themselves are 
sampled in a degree-corrected manner, conditioned on 
the edge labels. 


P{{id)\{r,s),K,c) 



, r ^ Cj , s 
2i ,T^Cj 


if r 7 ^ s 
if r = s. 


(14) 


where is the probability of a node being selected in¬ 
side a group, with total likelihood 

conditioned on the label sequence becomes 


use noninformative priors, but this time on {Ki}, inte¬ 
grate over them, obtaining (omitting henceforth the triv¬ 
ial Kronecker delta term above) 


TT fJ ^ TT 

i^({(bj)*}|{(D5)a,c) = -P({dJ), (16) 


n. 


with P{{di}) = J([,. (( • Combining this with Eq. 

as P{{{i,j)t}\c, b) = P({(l j)t}|{(r, s)t}, c)P({(r, s)t}\b)^ 
we have the complete likelihood of the temporal network 


Pi{ihj)t}\c,b) = 


nr>s"ir-s!nr- 2" 


Ur er-! 




xP{{di}\c)P{{mr.s}) 


n 


6^ I 

U<V ^UV' 


n« pu-Uv pv'- 


jlPi{eL})- (18) 


This likelihood can be rewritten in such a way that makes 
clear that it is composed of one purely static and one 
purely dynamic part. 


P{{{t,j)t}\c,b) = P{{A,}\c)x 


-P({(Dg)t}IMe.;}) 


(19) 


The first term of Eq. is precisely the nonparametric 
likelihood of the static DCSBM that generates the aggre¬ 
gated graph with adjacency matrix Aij = kx=(ij) given 
the node partition {q}, which itself is given by 


lnP({Ay}|c) « E -I- i ^ e^s In ^ IndJ 

^ Gr Gc 

rs i 

-I-lnP({di})-I-lnP({mrs}), (20) 


P{{ihj)t}\{{r,s)t},K,c) = ]JP((*, j)t|(r,s)t,K) 




(15) 


where di is the degree of node i, and m^s is the total 
number of edges between groups r and s. Performing 
maximum likelihood, one obtains ki = di/Cc,- But since 
we want to avoid overfitting the model, we once more 


if Stirling’s approximation is used. The second term in 
Eq.[^is the likelihood of the Markov chain of edge labels 
given by Eq.[^ (with {xt} = {(r, s)t}, and = {rrirs})- 
This model, therefore, is a direct generalization of the 
static DCSBM, with a likelihood composed of two sepa¬ 
rate static and dynamic terms. One recovers the static 
DCSBM exactly by choosing Bm = Bm = 1 — making 
the state transitions memoryless — so that the second 
term in Eq. above contributes only with a trivial con¬ 
stant 1/E\ to the overall likelihood. Equivalently, we can 
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Table II. Description length S = — in P({(i, c, 6) (in nats) as well as inferred number of node, token and memory groups, 
C, Bn and Bm, respectively, for different data sets and different Markov order n (see Appendix 0 - The value —AS < 
in P({a;(}|{a;*}, 6*) is a lower-bound on the predictive likelihood of the validation set {x'{\ (corresponding to half of the entire 
sequence) given the training set {x*} and its best parameter estimate. 



Figure 4. Inferred temporal model for a proximity network of 
high-school students m- (a) The static part of the model di¬ 
vides the students into C = 10 groups (square nodes) that al¬ 
most match the known classes (text labels). (b) The dynamic 
part of the model divides the directed multigraph group pairs 
in (a) into Bn = Bm = 9 groups (grey circles). The model 
corresponds to a n = 1 unified Markov chain on the edge la¬ 
bels, where the memory and tokens have identical partitions 
(see Appendix [ b| for details). 


view the DCSBM as a special case with n = 0 of this 
temporal network model. 

To make the model nonparametric, we again in¬ 
clude the same prior as before for the node partition 
c, in addition to token/memory partition b, such that 
the total nonparametric joint likelihood is maximized, 
= P{{{i,j)t}\c,b)P{c)Pib). In this way 
we are again protected against overfitting, and we can 
infer not only the number of memory and token groups, 
Bn and Bm, as before, but also the number of groups 
in the temporal network itself, C. If, for example, the 
temporal network is completely random — i.e. the edges 
are placed randomly both in the aggregated network as 
well as in time — we always infer Bn = Bm = C = 1. 

We employ this model in a variety of dynamic network 
datasets from different domains (see Table [IT] and Ap¬ 


pendix]^ for dataset descriptions). In all cases, we infer 
models with n > 0 that identify many groups for the 
tokens and memories, meaning that the model succeeds 
in capturing temporal structures. In most cases, mod¬ 
els with n = 1 best describe the data, implying there is 
no sufficient evidence for higher-order memory, with the 
exception of the network of chess moves, which is best 
described by a model with n = 2. To illustrate how the 
model characterizes the temporal structure of these sys¬ 
tems, we focus on the proximity network of high school 
students, which corresponds to the voluntary tracking of 
327 students for a period of 5 days m- Whenever the 
distance between two students fell below a threshold, an 
edge between them was recorded at that time. In Fig. 
we see the best-fitting model for this data. The groups 
inferred for the aggregated network correspond exactly 
to the known division into 9 classes, as indicated in the 
figure, with the exception of the PC class, which was 
divided into two groups. The groups show a clear as- 
sortative structure, where most connections occur within 
each class. The clustering of the edge labels in the sec¬ 
ond part of the model reveals the temporal dynamics. 
We observe that the edges connecting nodes of the same 
group cluster either in single-node or small groups, with 
a high incidence of self-loops. This means that if an edge 
that connects two students of the same class appears in 
the sequence, the next edge is most likely also inside the 
same class, indicating that the students of the same class 
are clustered in space and time. The remaining edges be¬ 
tween students of different classes are separated into two 
large groups. This division indicates that the different 
classes meet each other at different times. Indeed, the 
classes are located in different parts of the school build¬ 
ing and typically go to lunch separately m- 
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A. Temporal prediction 

A direct advantage of being able to extract such tem¬ 
poral patterns is that they can be used to make pre¬ 
dictions. This is in particular true of the Bayesian ap¬ 
proach, since it can even be used to predict tokens and 
memories not previously observed. We demonstrate this 
by dividing a sequence into two equal-sized contiguous 
parts, {xt} = {xl} U {x^}, i.e. a training set {a;^} and a 
validation set {x'^}. If we observe only the training set, 
a lower bound on likelihood of the validation set con¬ 
ditioned on it given by P{{x^}\{xl},b*) > exp(—AS) 
where b' = aTgmaXf^, P{{xf} L) {Xf}\b* ,b')P{b* ,b') and 
AS is the difference in the description length between 
the training set and the entire data (see Appendix 
for a proof). This lower bound will become tight when 
both the validation and training sets become large, and 
hence can be used as an asymptotic approximation of the 
predictive likelihood. Table [II] shows empirical values for 
the same datasets as considered before, where n = 0 cor¬ 
responds to using only the static DCSBM to predict the 
edges, ignoring any time structure. The temporal net¬ 
work model provides better prediction in all cases. 


IV. MODEL EXTENSIONS 
A. Continuous time 

So far, we have considered sequences and temporal net¬ 
works that evolve discretely in time. Although this is the 
appropriate description for many types of data, such as 
text, flight itineraries and chess moves, in many other 
cases events happen instead in real time. In this case, 
the time series can be represented — without any loss of 
generality — by an “embedded” sequence of tokens {xt} 
placed in discrete time, together with an additional se¬ 
quence of waiting times {A^}, where Aj > 0 is the real 
time difference between tokens xt and Xt-i- Employing a 
continuous-time Markov chain description, the data like¬ 
lihood can be written as 

P{{xt},{At}\p,X) = P{{xt}\p) X P({At}|{a:t}, A) (21) 
with P{{xt\\p) given by Eq. and 

P({Aa|K},A) (22) 


where 


P(A|A) = Ae"^'^, (23) 

is a maximum-entropy distribution governing the waiting 
times, according to the frequency A. Substituting this in 
Eq. we have 

P({AJ|{cra, A) = n (24) 


where A^j = To compute the nonparamet- 

ric Bayesian evidence, we need a conjugate prior for the 
frequencies Xg, 


P(A|a,/3) 


r(a) 


(25) 


where a and /3 are hyperparameters, interpreted, respec¬ 
tively, as the number and sum of prior observations. 
A fully noninformative choice would entail a —>■ 0 and 
/? —)■ 0, which would yield the so-called Jeffreys prior im, 
P(A) oc 1/A. Unfortunately, this prior is improper, i.e. it 
is not a normalized distribution. In order to avoid this, 
we use instead a = 1 and /3 = 'YhxXs/XI, taking into 
account the global average. Using this prior, we obtain 
the Bayesian evidence for the waiting times as 




poo 

n/ dAA'=^-e-"^--p(A|a,/3), 

X “'o 

-j-r /3°T(fcg -I- a) 

-lir(a)(As + / 3 )fc-+“' 


(26) 

(27) 


Hence, if we employ the Bayesian parametrization with 
communities for the discrete embedded model as we did 
previously, we have 


P{{xt}, {AJ, b) = P{{xt}, b) X P{{At}\{xt}), (28) 

with P{{xt},b) given by Eq. [I^ 

Since the partition of memories and tokens only influ¬ 
ences the first term of Eq. |28| corresponding to the em¬ 
bedded discrete-time Markov chain, P{{xt},b), the out¬ 
come of the inference for any particular Markov order will 
not take into account the distribution of waiting times — 
although the preferred Markov order might be influenced 
by it. We can change this by modifying the model above, 
assuming that the waiting times are conditioned on the 
group membership of the memories, 

X3 = Vbs, (29) 


where rjr is a frequency associated with memory group 
r. The Bayesian evidence is computed in the same man¬ 
ner, integrating over rjj. with the noninformative prior of 
Eq. yielding 


Pi{\}\{xt}) 


TT /3°T(er -I- a) 

-l-J- r(a)(Ar -t- /3)®-+“’ 


(30) 


where A^ = J2t ^t5bs^,v 

As an example of the use of this model variation, 
we consider a piano reduction of Beethoven’s fifth sym¬ 
phony |52], represented as a sequence of A = 4, 223 notes 
of an alphabet of size N = 63. We consider both model 
variants where the timings between notes are discarded, 
and where they are included. If individual notes occur 
simultaneously as part of a chord, we consider them to 
have occurred in relative alphabetical order, with a in¬ 
terval of At = 10“® seconds between them. The results 





Continuous time 


Discrete time 


n 

Bn 

Bm 

- lnP({a;t}, b) 

Bn 

Bm 

-lnP({a;t},{At},6) 

1 

40 

40 

13,736 T 

37 

37 

58,128 

2 

35 

34 

15,768 

24 

22 

47,408 1 

3 

34 

33 

24, 877 

16 

15 

54,937 


Table III. Joint likelihoods for the discrete- and continuous¬ 
time Markov models, inferred for a piano reduction of 
Beethoven’s fifth symphony, as described in the text, for dif¬ 
ferent Markov orders n. Values in grey correspond to the 
maximum likelihood of each column. 


of the inference can be seen in table cni The discrete¬ 
time model favors a n = 1 Markov chain, whereas the 
continuous-time model favors n = 2. This is an interest¬ 
ing result that shows that the timings alone can influence 
the most appropriate Markov order. We can see in more 
detail why by inspecting the typical waiting times condi¬ 
tioned on the memory groups, as shown in Fig. For the 
discrete-time model, the actual continuous waiting times 
(which are not used during inference) are only weakly 
correlated with the memory groups. On the other hand, 
for the continuous-time model we And that the memories 
are divided in such a way that they are very strongly 
correlated with the waiting times: There is a group of 
memories for which the ensuing waiting times are always 
At = 10“® — corresponding to node combinations that 
are always associated with chords. The remaining mem¬ 
ories are divided into further groups that display at least 
two distinct time scales, i.e. short and long pauses be¬ 
tween notes. These statistically significant patterns are 
only visible for the higher order n = 2 model. 


1. Bursty dynamics 
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Figure 5. Average waiting time per memory group, (A)^ = 
Ar/rir, for the Markov models inferred from Beethoven’s fifth 
symphony, (a) n = 1 discrete-time model, ignoring the wait¬ 
ing times between notes, (b) n = 2 continuous-time model, 
with waiting times incorporated into the inference. 


In the above model the waiting times are distributed 
according to the exponential distribution of Eq.|^ which 
has a typical time scale given by 1/A. However, one often 
encounters processes where the dynamics is bursty, i.e. 
the waiting times between events lack any characteristic 
scale, and are thus distributed according to a power-law 

P{A\B) = 1^, (31) 

for A > Am, otherwise P{A\B) = 0. One could in princi¬ 
ple repeat the above calculations with the above distribu¬ 
tion to obtain the inference procedure for this alternative 
model. However, this is in fact not necessary, since by 
making the transformation of variable 

p = ln^, (32) 

we obtain for Eq. 

P{fi\B) = Be^^, (33) 

which is the same exponential distribution of Eq. |23| 
Hence, we need only to perform the transformation of 


Eq.j^for the waiting times prior to inference, to use the 
bursty model variant, while maintaining the exact same 
algorithm. 


B. Nonstationarity 

An underlying assumption of the Markov model pro¬ 
posed is that the same transition probabilities are used 
for the whole duration of the sequence, i.e. the Markov 
chain is stationary. Generalizations of the model can be 
considered where these probabilities change over time. 
Perhaps the simplest generalization is to assume that the 
dynamics is divided into T discrete “epochs”, such that 
one replaces tokens Xt by a pair (x,T)t, where r S [1,T] 
represents the epoch where token x was observed. In 
fact, T does not need to be associated with a temporal 
variable — it could be any arbitrary covariate that de¬ 
scribes additional aspects of the data. By incorporating 
this type of “annotation” into the tokens, one can use 
a stationary Markov chain describing the augmented to¬ 
kens that in fact corresponds to a non-stationary one if 
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Tokens Memories Tokens Memories 


(a) 

Figure 6. Markov model fit for a concatenated text composed 
of “War and peace”, by Leo Tolstoy, and “A la recherche du 
temps perdu”, by Marcel Proust. Edge endpoints in red and 
blue correspond to token an memories, respectively, that in¬ 
volve letters exclusive to French, (a) Version of the model 
with n = 3 where no distinction in made between the same 
token occurring in the different novels, yielding a description 
length — log 2 P{{xt}, b) = 7,450,322. (b) Version with n = 3 
where each token is annotated with its occurring novel, yield¬ 
ing a description length — log 2 P{{xt}, b) = 7,146,465. 


one omits the variable r from the token descriptors — ef¬ 
fectively allowing for arbitrary extensions of the model by 
simply incorporating appropriate covariates, and without 
requiring any modification to the inference algorithm. 

Another consequence of this extension is that the same 
token X can belong to different groups if it is associated 
with two or more different covariates, {x,ti) and (x,T 2 ). 
Therefore, this inherently models a situation where the 
group membership of tokens and memory vary in time. 

As an illustration of this application of the model, 
we consider two literary texts: an English translation of 
“War and peace”, by Leo Tolstoy, and the French original 
of “A la recherche du temps perdu”, by Marcel Proust. 


First, we concatenate both novels together, treating it as 
a single text. If we fit our Markov model to it, we obtain 
the n = 3 model shown in Fig. [§i. In that figure, we have 
highlighted tokens and memories that involve letters that 
are exclusive to the French language, and thus occur most 
predominantly in the second novel. We observe that the 
model essentially finds a mixture between English and 
French. If, however, we indicate in each token to which 
novel it belongs, e.g. (a:,wp)t and (a;,temps)i, we obtain 
the model of Fig. In this case, the model is forced to 
separate between the two novels, and one indeed learns 
the French patterns differently from English. Since this 
“nonstationary” model possesses a larger number of mem¬ 
ory and tokens, one would expect a larger description 
length. However, in this cases it has a smaller descrip¬ 
tion length than the mixed alternative, indicating indeed 
that both patterns are sufficiently different to warrant a 
separate description. Therefore, this approach is capable 
of uncovering change points US], where the rules govern¬ 
ing the dynamics change significantly from one period to 
another. 


V. CONCLUSION 

We presented a dynamical variation of the degree- 
corrected stochastic block model that can capture long 
pathways or large-scale structures in sequences and tem¬ 
poral networks. The model is based on a nonparametric 
variable-order hidden Markov chain, and can be used not 
only to infer the most appropriate Markov order but also 
the number of groups in the model. Because of its non¬ 
parametric nature, it requires no a priori imposition of 
time scales, which get inferred according to the statistical 
evidence available. 

We showed that the model successfully finds meaning¬ 
ful large-scale temporal structures in empirical systems, 
and that it predicts their temporal evolution. We demon¬ 
strated also the versatility of the approach, by showing 
how the model can be easily extended to situations where 
the dynamics occur in continuous time, or is nonstation¬ 
ary. 

Our approach provides a principled and more power¬ 
ful alternative to approaches that force the division of 
the network-formation dynamics into discrete time win¬ 
dows, and require the appropriate amount of dynamical 
memory to be known in advance. 
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Appendix A: Bayesian Markov chains with 
commnnities 


As described in the main text, a Bayesian formulation 
of the Markov model consists in specifying prior prob¬ 
abilities for the model parameters, and integrating over 
them. In doing so, we convert the problem from one 
of parametric inference (i.e. when the model parameters 
need to be specified or fitted from data), to a nonpara- 
metric one (i.e. no parameters need to be specified or 
fitted). In this way, the approach possesses intrinsic reg¬ 
ularization, where the order of the model can be inferred 
from data alone, without overfitting iSHUsg. 

To accomplish this, we rewrite the model likelihood 
(using Eqs. and as 

Pi{xt}\b,x,9 )> 

x,x X r<s 

(Al) 

and observe the normalization constraints Y) - 0^ = 1, 
and Xrs = 1- Since this is just a product of multino¬ 
mials, we can choose conjugate Dirichlet priors probabil¬ 
ity densities 'Dr{{9x}\{ax}) and X>s({Ars}|{/3r j), which 
allows us to exactly compute the integrated likelihood. 


P{{xt}\a,l3,b) = J dOdX Pi{xt}\b,X,9) 


r s 

n r(Ar) -pr r(A;a; J-Qfj;) 

Tier + Ar) y r(a,) 

\_r xGr _ 

r(iJs) -r-j- r(ers +/3rs) 

^ lir(e, + i?,)li r(/3,.,) 


(A2) 
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where Ar = J^xer = Sr ^'rs- We recover the 

Bayesian version of the common Markov chain formu¬ 
lation (see Ref. [33]) if we put each memory and token 
in their own groups. This remains a parametric distri¬ 
bution, since we need to specify or infer the hyperpa¬ 
rameters. However, in the absence of prior information 
it is more appropriate to make a noninformative choice 
that encodes our a priori lack of knowledge or preference 
towards any particular model, (Xx = firs = 1 , which sim¬ 
plifies the above in a way that can be written exactly 
as Eqs. 6-9 in the main text. As was mentioned there, 
the likelihood of Eq. 7 corresponds to a microcanonical 
model that generates random sequences with hard con¬ 
straints. In order to see this, consider a chain where there 
are only e^s transitions in total between token group r 
and memory group s, and each token x occurs exactly kx 
times. For the very first transition in the chain, from a 
memory xq in group s to a token xi in group r, we have 
the probability 

P{xi\xo,b) = . (A3) 

Now, for the second transition from memory xi in group 
t to a token X 2 in group u, we have the probability 

P{x2\xi,b) = 

if t yf s, u yf r, X 2 ^ Xi, 

if t = s, u ^ r, X 2 ^ Xi, 

if t yf s, u = r, X 2 = Xi, 

if t yf s, u = r, X 2 ^ Xi, 

if t = s, u = r, X 2 ^ Xi, 

if t = s, u = r, X 2 = Xi- 

(A4) 

Proceeding recursively, the final likelihood for the entire 
chain is 

P{{Xt}\bAers},{kx}) = n—fw" , (A5) 

llr er irises! Y 

which is identical to Eq. 7 in the main text. 

Since the integrated likelihood above gives 
P{{xt}\b^{es}), we still need to include priors for the 
node partitions {bx} and {b^}, as well as memory group 
counts, {bs}, to make the above model fully nonpara- 
metric. This is exactly the same situation encountered 
with the SBM |31| |13J [H]. Following Refs. (33] lll|, we 
use a nonparametric two-level Bayesian hierarchy for 
the partitions, P{{bi}) = P{{bi}\{nr})P{{nr}), with 
uniform distributions 

^(MIK}) = S^, P{{nr})={^Z^^ \ 

(A6) 


^ut^X2 

\^US ^)^X2 
1) 

eticr - 1 ) 

^rtkx2 
Bt(Br 1) 

(Ors f^kx2 
es - l)(e,. - 1)’ 
Brs f^{kxi 1) 
. (Bs - l){er - 1) 


where M = tir, which we use for both {b^} and {b^}, 
i.e. P{b) = P{{bx})P{{bs;}). Analogously, for {bs} we 
can use a uniform distribution 

P{{es}\b)=(I^^^Jj \ (A7) 

The above priors make the model fully nonparamet¬ 
ric with a joint/marginal probability P{{xt},b) = 

P{{xt},b,{es}) = P{{xt}\b,{ee})P{b)P{{es}). As has 
been shown in Refs. |43fl33| . this approach succeeds in 
finding the most appropriate model size (i.e. number of 
groups) according to statistical evidence, without over¬ 
fitting. And also like in Ref. |33|, this nonparamet¬ 
ric method can be used to detect the most appropri¬ 
ate order of the Markov chain, again without overfit- 
ting. However, in some ways it is still sub-optimal. The 
use of conjugate Dirichlet priors above was primarily for 
mathematical convenience, not because they closely rep¬ 
resent the actual mechanisms believed to generate the 
data. Although the noninformative choice of the Dirich¬ 
let distribution (which yields flat priors for flat priors 
for {crs} and {bs}) can be well justified via maximum 
entropy arguments (see Ref. M), and are unbiased, it 
can in fact be shown that it can lead to underfitting 
of the data, where the maximum number of detectable 
groups scales sub-optimally as VJV [33]. As shown in 
Ref. |44|, this limitation can be overcome by depart¬ 
ing from the model with Dirichlet priors, and replac¬ 
ing directly the priors P({ei.s}|{es}) and P({ej;}) of the 
microcanonical model by a single prior P({ers}), and 
noticing that {e^s} corresponds to the adjacency ma¬ 
trix of bipartite multigraph with P edges and Bjv + Bm 
nodes. With this insight, we can write P{{ers}) as a 
Bayesian hierarchy of nested SBMs, which replaces the 
resolution limit above by N/ In N, and provides a multi¬ 
level description of the data, while remaining unbiased. 
Furthermore, the uniform prior in Eq. 8 for the token 
frequencies P{{kx}\{ers},b) intrinsically favors concen¬ 
trated distributions of values. Very often (e.g. in 
text and networks) this distribution is highly skewed. 
We therefore replace it by a two-level Bayesian hierarchy 
P{{kx}\{ers}, b) = rir Pi{kx}\{nl})Pi{nl}\er), with 

P{{kx}\{nl}) = (A8) 

and Pi{n^}\er) = q{er,nr)~^, where q{m,n) is the num¬ 
ber of restricted partitions of integer m into at most n 
parts (see Ref. [39] for details). 

As mentioned in the main text, in order to fit the model 
above we need to find the partitions {b^} and {bg} that 
maximize P{{xt},b), or fully equivalently, minimize the 
description length E = — lnP({a;t},6) [47]. Since this is 
functionally equivalent to inferring the DCSBM in net¬ 
works, we can use the same algorithms. In this work we 
employed the fast multilevel MCMC method of Ref. [3S] , 
which has log-linear complexity 0{N log^ N), where N is 
the number of nodes (in our case, memories and tokens), 
independent of the number of groups. 
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Appendix B: The unified n = 1 alternative model 

The model defined in the main text is based on a co¬ 
clustering of memory and tokens. In the n = 1 case, each 
memory corresponds to a single token. In this situation, 
we consider a slight variation of the model where we force 
the number of groups of each type to be the same, i.e. 
Bjy = Bm = B, and both partitions to be identical. 
Instead of clustering the original bipartite graph, this 
is analogous to clustering its projection into a directed 
transition graph with each node representing a specific 
memory and token simultaneously. When considering 
this model, the likelihoods computed in the main text 
and above remain exactly the same, with the only dif¬ 
ference that we implicitly force both memory and token 
partitions to be identical, and omit the partition likeli¬ 
hood of Eq. |A6| for one of them. We find that for many 
datasets this variation provides a slightly better descrip¬ 
tion than the co-clustering version, although there are 
also exceptions to this. 

We used this variation of the model in Fig. [^because 
it yielded a smaller description length for that dataset, 
and made easier the visualization and interpretation of 
the results in that particular case. 


Appendix C: Predictive held-out likelihood 


Given a sequence divided in two contiguous parts, 
{xt} = {xl} U {x[}, i.e. a training set {a;^} and a val¬ 
idation set {x[}, and if we observe only the training set, 
the predictive likelihood of the validation set conditioned 
on it is 


p{{x'}\{xn,b*) 


P{{x',}U{xt}\b*) 

p{{xn\b*) 


(Cl) 


where b* = argmax^ P(6|{a;t}) is the best partition given 
the training set. In the above, we have 


P{{x',} U {xtW) = ^ P{{x't} U {xn\b*,b')P{b'\b*), 

(C2) 

where b' corresponds to the partition of the newly ob¬ 
served memories (or even tokens) in {a;(}. Generally we 
have P{b'\b*) = P{b',b*)/P{b*), so that 


P({x;}|K},n = 


Y.,,P{{x',}VJ{xtW,b')P{P,b') 


> 


P{{xt}\b*)P{b*) 

p{{x',}^{xt}\b\h')p{b\y) 


P{{xl}\b*)P{b*) 


= exp(—AS), (C3) 


where b' = argmax^ P({a:(} U 6')P(&*, 6') and 

AE is the difference in the description length between 
the training set and the entire data. Hence, computing 
the minimum description length of the remaining data, 
via a maximization of the posterior likelihood relative to 
the partition of the previously unobserved memories or 


tokens, yields a lower bound on the predictive likelihood. 
This lower bound will become tight when both the vali¬ 
dation and training sets become large, because then the 
posterior distributions concentrate around the maximum, 
and hence can be used as an asymptotic approximation 
of the predictive likelihood. 


Appendix D: Equivalence between structure and 
dynamics 

The likelihood of Eq. 4 in the main text is almost the 
same as the DCSBM [3S]. The only exceptions are triv¬ 
ial additive and multiplicative constants, as well as the 
fact that the degrees of the memories do not appear in 
it. These differences, however, do not alter the position 
of the maximum with the respect to the node partition. 
This allows us to establish an equivalence between in¬ 
ferring the community structure of networks and mod¬ 
elling the dynamics taking place on it. Namely, for a 
random walk on a connected undirected graph, a tran¬ 
sition j —>■ j is observed with probability AijPi(t)/ki, 
with Pi(t) being the occupation probability of node i at 
time t. Thus, after equilibration with Pi{oo) = ki/2E, 
the probability of observing any edge (i,j) is a constant: 
Pi{oo)/ki +pj{oo)/kj = 1/E. Hence, the expected edge 
counts Crs between two groups in the Markov chain will 
be proportional to the actual edge counts in the under¬ 
lying graph given the same node partition. This means 
that the likelihood of Eq. 4 in the main text (for the n = 1 
projected model described above) and of the DGSBM will 
differ only in trivial multiplicative and additive constants, 
such that the node partition that maximizes them will be 
identical. This is similar to the equivalence between net¬ 
work modularity and random walks |22| . but here the 
equivalence is stronger and we are not constrained to 
purely assortative modules. However, this equivalence 
breaks down for directed graphs, higher order memory 
with n > 1 and when model selection is performed to 
choose the number of groups. 


Appendix E: Comparison with the map equation for 
network fiows with memory 

Both the community-based Markov model introduced 
here and the map equation for network fiows with mem¬ 
ory |24| identify communities in higher-order Markov 
chains based on maximum compression. However, the 
two approaches differ from each other in some central 
aspects. The approach presented here is based on the 
Bayesian formulation of a generative model, whereas the 
map equation finds a minimal entropy encoding of the 
observed dynamics projected on a node partition. Thus, 
both approaches seek compression, but of different as¬ 
pects of the data. 

The map equation operates on the internal and ex¬ 
ternal transitions within and between possibly nested 
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groups of memory states and describes the transitions be¬ 
tween physical nodes \xt is the physical node or token in 
memory states of the form x = {xt,Xt-i,Xt- 2 , ■ ■ ■)]• The 
description length of these transitions is minimized for 
the optimal division of the network into communities. By 
constrnction, this approach identifies assortative modnles 
of memory states with long flow persistence times. More¬ 
over, for inferring the most appropriate Markov order, 
this dynamics approach reqnires supervised approaches 
to model selection that uses random subsets of the data 
such as bootstrapping or cross validation |54| . 

On the other hand, the model presented here yields a 
nonparametric log-likelihood for the entire sequence as 
well as the model parameters, with its negative valne 
corresponding to a description length for the entire data, 
not only its projection into groups. The minimization of 
this description length yields the optimal co-clustering of 
memories and tokens, and hence no inherent assortativ- 
ity is assumed. Therefore it can be nsed also when the 
nnderlying Markov chain is dissortative. The description 
length also can be nsed to perform nnsnpervised model 
selection, where the Markov order and number of groups 
are determined from the entire data, obviating the need 
for bootstrapping or cross validation. Furthermore, since 
after the inference we have a trained generative model, 
the present approach can be used to generate new data 
and make predictions, based on past observations. 

Becanse of the above distinctions, the two different ap¬ 
proaches can give different resnlts and the problem at 
hand shonld decide which method to use. 

Appendix F: Datasets 

Below we give a description of the datasets used in this 
work. 

US Air Flights: This dataset corresponds to a sample 
of air flight itineraries in the US dnring 2011 col¬ 
lected by Bnreau of Transportation Statistics of the 
United States Department of Transportation |55| . 
The dataset contains 1, 272,696 itineraries of varied 
lengths (airport stops). We aggregate all itineraries 
into a single seqnence by concatenating the individ- 
nal itineraries with a special separator token that 
marks the end of a single itinerary. There are 464 
airports in the dataset, and hence we have an al¬ 
phabet of N = 465 tokens, and a single seqnence 
with a total length of 83, 653,994 tokens. 

War and Peace: This dataset corresponds to the entire 
text of the english translation of the novel War and 
Peace by Leo Tolstoy, made available by the Project 
Gutenberg [^. This corresponds to a seqnence 
with an alphabet of size A = 84 (including letters, 
space, punctuation and special symbols) and a total 
length of 3, 226,652 tokens. 

Taxi movements: This dataset contains GPS logs from 
25,000 taxi pickups in San Francisco, collected by 


the company Uber m- The geographical loca¬ 
tions were discretized into 416 hexagonal cells (see 
Ref. [23] for details), and the taxi rides were con¬ 
catenated together in a single sequence with a spe¬ 
cial separator token indicating the termination of 
a ride. In total, the sequence has an alphabet of 
N = 417 and a length of 819,172 tokens. 

RockYou password list: This dataset corresponds to 
a widely distributed list of 32, 603, 388 passwords 
from the RockYou video game company |5S|. The 
passwords were concatenated in a single sequence, 
with a special separator token between passwords. 
This yields a sequence with an alphabet of size 
N = 215 (letters, numbers and symbols) and a to¬ 
tal length of 289,836, 299 tokens. 

High school proximity: This dataset corresponds to a 
temporal proximity measnrement of stndents in a 
french high school m- A total of A = 327 students 
were voluntarily tracked for a period of five days, 
generating if = 5,818 proximity events between 
pairs of stndents. 

Enron email: This dataset corresponds to time- 
stamped collection of if = 1,148,072 emails 

between directed pairs of A = 87,273, senders 
and recipients of the former Enron corporation, 
disclosed as part of a frand investigation |59| . 

Internet AS: This dataset contains connections be¬ 
tween autonomons systems (AS) collected by the 
GAIDA project m- It corresponds to a time- 
stamped seqnence of if = 500,106 directed connec¬ 
tions between AS pairs, with a total of A = 53, 387 
recorded AS nodes. The time-stamps correspond 
to the first time the connection was seen. 

APS citations: This dataset contains E = 4,262,443 
time-stamped citations between A = 425, 760 sci¬ 
entific articles published by the American Physical 
Society for a period of over 100 years [6T| 

prosper.com loans: This dataset corresponds to if = 
3,394,979 time-stamped directed loan relation¬ 
ships between A = 89,269 nsers of the prosper, 
com website, which provides a peer-to-peer lending 
system. |H2] 

Chess moves: This dataset contains 38,365 online 
chess games collected over the month of July 
2015 [B3|- The games were converted into a bi¬ 
partite temporal network where each piece and po¬ 
sition correspond to different nodes, and a move¬ 
ment in the game corresponds to a time-stamped 
edge of the type piece —>■ position. The resnlting 
temporal network consists of A = 76 nodes and 
if = 3,130,166 edges. 

Hospital contacts: This dataset corresponds to a tem¬ 
poral proximity measurement of patients and 
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health care workers in the geriatric unit of an uni¬ 
versity hospital in]- A total of = 75 individuals 
were voluntarily tracked for a period of four days, 
generating E = 32,424 proximity events between 
pairs of individuals. 

Infectious Sociopatterns: This dataset corresponds 
to a temporal proximity measurement of attendants 
at a museum exhibition [5^. A total of iV = 10,972 
participants were voluntarily tracked for a period of 


three months, generating E = 415,912 proximity 
events between pairs of individuals. 

Reality Mining: This dataset corresponds to a tempo¬ 
ral proximity measurement of university students 
and faculty [5S]. A total of = 96 people were 
voluntarily tracked for a period of an entire aca¬ 
demic year, generating E = 1,086,404 proximity 
events between pairs of individuals. 


